#!/usr/bin/perl -w

###############################################################################
# Transform the annotation input files from the illumina platform
# into a database load file for the anno table in the database
#
# Parameters needed: dbname, username, anno_files, temp_dir
# Output files: temp_dir/anno.out
#
# Author: Herodotos Herodotou and Nedyalko Borisov
# Date: March 11, 2010
#
##############################################################################

use getOpts qw(get_options println);
use allele_convertion qw(convert_alleles flip_allele get_chromosome_number);

if ($#ARGV != 0)
{
   println qq(Please specify the .properties file);
   println qq(Usage: perl $0 load.properties);
   exit(-1);
}

# Get input parameters
%opt = get_options($ARGV[0]);
$DB_NAME  =$opt{'dbname'};
$USERNAME =$opt{'username'};
$ANNO_FILES = $opt{'anno_files'};
$TEMP_DIR = $opt{'temp_dir'};

# Open output file
open OUTFILE, ">", $TEMP_DIR."/anno.out" or die $!;

# Create the snp reference file
$query = qq(SELECT rsid FROM anno WHERE platform = 'i');
!system qq(psql -U $USERNAME -d $DB_NAME ).
        qq(-c "$query" -t -o $TEMP_DIR/anno_ref) or die $!;

# Load the reference file in the anno map
open INFILE, "<", "$TEMP_DIR/anno_ref" or die $!;

%snps = ();
while ($line = <INFILE>)
{
   $line =~ s/^\s+//;
   $line =~ s/\s+$//;
   $snps{$line} = '';
}

close INFILE;
!system qq(rm -f $TEMP_DIR/anno_ref) or die $!;


# Perform the transformation
@annofiles = split(/;/, $ANNO_FILES);

foreach (@annofiles)
{
    $annofile = $_;
    $annofile =~ s/\s+$//;
    $num_records = 0;

    # Open input anno files
    println qq(Transform anno file: $annofile);
    open INFILE, "<", $annofile or die $!;

    # Find the header line
    while($line = <INFILE>)
    {
       if ($line =~ "^IlmnID")
       {
          last;
       }
    }

    # Process the data
    while($line = <INFILE>)
    {
       # Check for the end of the actual data
       if ($line =~ "^\\[Controls\\]")
       {
          last;
       }

       # Get the data
       @data = split(/,/, $line);
       $manuf_id = qq($data[0]);
       $rsid = $data[1];
       $chromosome = get_chromosome_number($data[9]);
       $location = ($data[10] =~ /---/) ? "-9" : $data[10];

       # Make sure this snp is not already in the database   
       if (exists $snps{$rsid})
       {
          next;
       }

       # Make the necessary allele convertions
       $alleleA = substr($data[3], 1, 1);
       $alleleB = substr($data[3], 3, 1);   
       ($alleleA, $alleleB, $reverse) 
           = convert_alleles($data[2], $data[15], $alleleA, $alleleB);
           
       if (!defined($alleleA) || !defined($alleleB) || !defined($reverse))
       {
          println qq(Incorrect format for line: $line);
          next;
       }
       
       # Print anno
       print OUTFILE qq($manuf_id,i,$rsid,$chromosome,$location,$alleleA,$alleleB,$reverse\n);
       ++$num_records;
    }
    
    # Close file
    close INFILE;
    println qq(Done anno file transformation - produced $num_records records);
}

# Done
close OUTFILE;

